Spatial Wrangling

GESIS Workshop: Introduction to Geospatial Techniques for Social Scientists in R

Stefan Jünger & Dennis Abel

2025-04-09

Now

Day Time Title
April 09 10:00-11:30 Introduction
April 09 11:30-11:45 Coffee Break
April 09 11:45-13:00 Data Formats
April 09 13:00-14:00 Lunch Break
April 09 14:00-15:30 Mapping I
April 09 15:30-15:45 Coffee Break
April 09 15:45-17:00 Spatial Wrangling
April 10 09:00-10:30 Mapping II
April 10 10:30-10:45 Coffee Break
April 10 10:45-12:00 Applied Spatial Linking
April 10 12:00-13:00 Lunch Break
April 10 13:00-14:30 Spatial Autocorrelation
April 10 14:30-14:45 Coffee Break
April 10 14:45-16:00 Spatial Econometrics & Outlook

Mini wrap-up

Thus far, we have learned about

  • different data formats
  • how to load them
  • first steps in interacting with them
  • creating maps with tmap

But the real magic in geospatial data lies in their flexibility.

Our plan for this afternoon

In this session, you will learn

  • how to wrangle the different geospatial data formats even further
  • how to converse from one format into the other

Because we want to

  • link/join different datasets
  • do some spatial analysis

1. Advanced Importing

Importing of non-spatial data

Say, the data we want to use are not available as a shapefile. Point coordinates are often stored in table formats like .csv – as is the location of charging stations for electric cars data in our ./data folder.

echarging_df <- 
  readr::read_delim("./data/charging_points_ger.csv", delim = ";")

head(echarging_df)
# A tibble: 6 × 7
  operator             federal_state latitude longitude power_kw type  num_plugs
  <chr>                <chr>            <dbl>     <dbl>    <dbl> <chr>     <dbl>
1 deer GmbH            Baden-Württe…     48.3      9.72       44 Norm…         2
2 EnBW mobility+ AG u… Baden-Württe…     48.6      9.87       93 Schn…         2
3 SWU Energie GmbH     Baden-Württe…     48.5     10.2        44 Norm…         2
4 SWU Energie GmbH     Baden-Württe…     48.6     10.1        44 Norm…         2
5 SWU Energie GmbH     Baden-Württe…     48.2     10.1        22 Norm…         1
6 EnBW mobility+ AG u… Baden-Württe…     48.5      9.98       30 Norm…         2

From data table to geospatial data

We see that besides our attributes (e.g., operator, power,…), the table contains the two variables “longitude” (X) and “latitude” (Y), our point coordinates. When using the command sf::st_as_sf(), it is easy to transform the table into a point layer.

# transform to spatial data frame
echarging_sf <- 
  sf::st_as_sf(
    echarging_df |>  
      # there were some missings in my data which is not allowed 
      dplyr::filter(!is.na(longitude) & !is.na(latitude)),    
    coords = c("longitude", "latitude")
  )

# inspect data
class(echarging_sf)
[1] "sf"         "tbl_df"     "tbl"        "data.frame"

Final check

plot(echarging_sf)

Check the CRS!

Make sure to use the option crs = [EPSG_ID]. If not used, your CRS will not be defined, and you can’t perform further commands depending on the CRS. Here, I tried EPSG IO or http://projfinder.com/ to find out.

echarging_sf <- 
  echarging_df |> 
  # there were some missings in my data which is not allowed  
  dplyr::filter(!is.na(longitude) & !is.na(latitude)) |> 
  sf::st_as_sf(   
    coords = c("longitude", "latitude"),
    crs = 4326
  )

plot(echarging_sf["operator"])

… and the other way round

Do you want to go back to handling a simple data frame? You can quickly achieve this by dropping the geometry column.

# remove geometry column
sf::st_drop_geometry(echarging_sf) |>  
  head(2)
# A tibble: 2 × 5
  operator                    federal_state     power_kw type          num_plugs
  <chr>                       <chr>                <dbl> <chr>             <dbl>
1 deer GmbH                   Baden-Württemberg       44 Normalladeei…         2
2 EnBW mobility+ AG und Co.KG Baden-Württemberg       93 Schnellladee…         2

APIs

Geospatial data tend to be quite big, there’s a pressure to distribute data efficiently. Data dumps (on the internet) may not be helpful

  • When resources are low
  • Time’s a factor
  • The data have a large geographic extent

Instead, a Programming Application Interface (API) is often used.

Data providers offering geospatial data APIs

Example: access to public transport

Say, we’re interested in the accessibility of public transport in Cologne

  • Bus, tram, etc.
  • All platforms and vehicles should be wheel-chair accessible

We can gather this information using OpenStreetMap!

Accessing OSM data: the Overpass API

The Overpass API (formerly known as OSM Server Side Scripting, or OSM3S before 2011) is a read-only API that serves up custom selected parts of the OSM map data. It acts as a database over the web: the client sends a query to the API and returns the data set that corresponds to the query.

Source: https://wiki.openstreetmap.org/wiki/Overpass_API

Starting with a geographic area to query

Many geospatial API requests start with a bounding box or another geographical extent to define which area should be accessed.

cologne_pt_stops <-
  osmdata::getbb(
    "Köln"
  )

cologne_pt_stops
       min       max
x  6.77253  7.162028
y 50.83044 51.084974

Defining some technical details

The Overpass API also requires a timeout parameter that repeats the request for a while if anything fails.

cologne_pt_stops <-
  cologne_pt_stops |> 
  osmdata::opq(timeout = 25*100)

cologne_pt_stops
$bbox
[1] "50.8304399,6.7725303,51.0849743,7.162028"

$prefix
[1] "[out:xml][timeout:2500];\n(\n"

$suffix
[1] ");\n(._;>;);\nout body;"

$features
NULL

$osm_types
[1] "node"     "way"      "relation"

attr(,"class")
[1] "list"           "overpass_query"
attr(,"nodes_only")
[1] FALSE

Turning to the content

The content we aim to request is defined with key/value pairs. It’s best to learn them by looking them up in the official documentation.

cologne_pt_stops <-
  cologne_pt_stops |>    
  osmdata::add_osm_feature(key = "public_transport", value = "stop_position")

cologne_pt_stops
$bbox
[1] "50.8304399,6.7725303,51.0849743,7.162028"

$prefix
[1] "[out:xml][timeout:2500];\n(\n"

$suffix
[1] ");\n(._;>;);\nout body;"

$features
[1] "[\"public_transport\"=\"stop_position\"]"

$osm_types
[1] "node"     "way"      "relation"

attr(,"class")
[1] "list"           "overpass_query"
attr(,"nodes_only")
[1] FALSE

Conduct actual request/download

The data is finally downloaded in the osmdata package, e.g., using the osmdata::osmdata_sf() function.

cologne_pt_stops <-
  cologne_pt_stops |>   
  osmdata::osmdata_sf()

cologne_pt_stops

Filter and transform

The data comprises a list that can be accessed as any list in R. Here, we extract the points and wrangle them (spatially).

cologne_pt_stops <-
  cologne_pt_stops$osm_points |> 
  tibble::as_tibble() |> 
  sf::st_as_sf() |> 
  dplyr::filter(wheelchair == "yes")

cologne_pt_stops
Simple feature collection with 593 features and 95 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 6.776977 ymin: 50.83381 xmax: 7.161797 ymax: 51.08398
Geodetic CRS:  WGS 84
# A tibble: 593 × 96
   osm_id   name    FIXME `VRS:gemeinde` `VRS:name` `VRS:old_ref` `VRS:ortsteil`
 * <chr>    <chr>   <chr> <chr>          <chr>      <chr>         <chr>         
 1 361716   Eifelp… <NA>  KÖLN           <NA>       <NA>          Innenstadt    
 2 388550   Brahms… <NA>  KÖLN           <NA>       <NA>          Lindenthal    
 3 21033643 Weinsb… <NA>  KÖLN           Weinsberg… <NA>          Ehrenfeld     
 4 27042401 Markt   <NA>  BERGISCH GLAD… Bergisch … <NA>          Mitte         
 5 28122005 Heumar… <NA>  KÖLN           <NA>       <NA>          Innenstadt    
 6 28301370 Kalker… <NA>  KÖLN           Kalker Fr… <NA>          Merheim       
 7 28301373 Merheim <NA>  KÖLN           <NA>       <NA>          Merheim       
 8 28301416 Refrath <NA>  BERGISCH GLAD… <NA>       <NA>          Refrath       
 9 30145575 Köln-H… <NA>  KÖLN           Holweide … <NA>          Holweide      
10 30388742 Mülhei… <NA>  <NA>           <NA>       <NA>          <NA>          
# ℹ 583 more rows
# ℹ 89 more variables: `VRS:ref` <chr>, alt_name <chr>, ashtray <chr>,
#   bench <chr>, bin <chr>, bus <chr>, bus_bay <chr>, bus_lines <chr>,
#   check_date <chr>, `check_date:crossing` <chr>, construction <chr>,
#   `construction:railway` <chr>, `construction:tram` <chr>, covered <chr>,
#   crossing <chr>, departures_board <chr>, description <chr>, direction <chr>,
#   disused <chr>, exit <chr>, fixme <chr>, `gtfs:name` <chr>, …

The data indeed are mappable

tm_shape(cologne_pt_stops) +
  tm_dots()

Example: GeoJSON files

JSON files are pretty popular

  • standardized and well-structured
  • similar to XML-files, but human readability is a bit better
  • also easy to parse for editors and browser

There’s also a JSON flavor for geospatial data…

GeoJSON snippet

{
"type": "FeatureCollection",
"features": [
{
"type": "Feature",
"id": 12,
"geometry": {
"type": "Polygon",
"coordinates": [
[
[
6.957362270020273,
50.94308762750329
]   
...

Source: https://www.offenedaten-koeln.de/

An application from Cologne’s Open Data portal

park_and_ride <-
  glue::glue(
    "https://geoportal.stadt-koeln.de/arcgis/rest/services/verkehr/",
    "verkehrskalender/MapServer/8/query?where=objectid+is+not+null&text=&",
    "objectIds=&time=&geometry=&geometryType=esriGeometryEnvelope&inSR=&",
    "spatialRel=esriSpatialRelIntersects&distance=&units=esriSRUnit_Foot&",
    "relationParam=&outFields=*&returnGeometry=true&returnTrueCurves=false&",
    "maxAllowableOffset=&geometryPrecision=&outSR=4326&havingClause=&",
    "returnIdsOnly=false&returnCountOnly=false&orderByFields=&",
    "groupByFieldsForStatistics=&outStatistics=&returnZ=false&returnM=false&",
    "gdbVersion=&historicMoment=&returnDistinctValues=false&resultOffset=&",
    "resultRecordCount=&returnExtentOnly=false&datumTransformation=&",
    "parameterValues=&rangeValues=&quantizationParameters=&featureEncoding=",
    "esriDefault&f=pjson"
  ) |> 
  sf::st_read(as_tibble = TRUE)
Reading layer `ESRIJSON' from data source 
  `https://geoportal.stadt-koeln.de/arcgis/rest/services/verkehr/verkehrskalender/MapServer/8/query?where=objectid+is+not+null&text=&objectIds=&time=&geometry=&geometryType=esriGeometryEnvelope&inSR=&spatialRel=esriSpatialRelIntersects&distance=&units=esriSRUnit_Foot&relationParam=&outFields=*&returnGeometry=true&returnTrueCurves=false&maxAllowableOffset=&geometryPrecision=&outSR=4326&havingClause=&returnIdsOnly=false&returnCountOnly=false&orderByFields=&groupByFieldsForStatistics=&outStatistics=&returnZ=false&returnM=false&gdbVersion=&historicMoment=&returnDistinctValues=false&resultOffset=&resultRecordCount=&returnExtentOnly=false&datumTransformation=&parameterValues=&rangeValues=&quantizationParameters=&featureEncoding=esriDefault&f=pjson' 
  using driver `ESRIJSON'
Simple feature collection with 24 features and 9 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 6.814912 ymin: 50.84826 xmax: 7.16169 ymax: 51.05211
Geodetic CRS:  WGS 84

Source: https://www.offenedaten-koeln.de/dataset/park-and-ride-anlagen-koeln

Park & ride parking spaces

With just two ‘simple’ commands, we can retrieve geospatial data about Cologne’s Park & Ride parking spaces in R. Not too bad, right?

tm_shape(park_and_ride) +
  tm_dots()

Raster data access

It is not only vector data that can be processed through these mechanisms.

The idea is the same for raster data

  • Accessing the information through URLs
  • Just the downloaded data formats differ

Example: downloading German weather data

The German Weather Service provides excellent weather data in its Climate Data Center (CDC): https://opendata.dwd.de/climate_environment/CDC/. Let’s download some summer temperature data using a custom function.

download_dwd_data <- function(url, path) {
  
  download.file(url, dest = paste0(path, "/", "lyr.1.asc.gz"))
  
  R.utils::gunzip(
    paste0(path, "/", "lyr.1.asc.gz"), 
    overwrite = TRUE,
    remove = TRUE
  )
  
  raster_file <-
    terra::rast(paste0(path, "/", "lyr.1.asc"))
  
  unlink(paste0(path, "/", "lyr.1.asc.gz"))
  
  raster_file
}

Voilà

paste0(
  "https://opendata.dwd.de/climate_environment/CDC/grids_germany/seasonal/",
  "air_temperature_max/14_JJA/grids_germany_seasonal_air_temp_max_202314.",
  "asc.gz") |> 
  download_dwd_data(path = "./data/") |> 
  terra::plot()

2. Linking

‘Spreadsheet join’

While a lot of our previous data were points, we only had a column for the German federal states as administrative information so far. We’re now moving “a layer down” and looking at Germany on a more fine-grained spatial level: the district. We just repeat what you already have done in the exercise 1_3_1_Basic Maps: joining data just like simple spreadsheets.

# load district shapefile
german_districts <- sf::read_sf("./data/VG250_KRS.shp")

# load district attributes
attributes_districts <- readr::read_csv2("./data/attributes_districts.csv") 

# join data
german_districts_enhanced <- 
  german_districts |>  
  dplyr::left_join(attributes_districts, by = "AGS")

Spatial join

But what can we do, if we do not have matching identifier. For example, in the German districts data and the e-charger data, there are no matching administrative identifier.

german_districts_enhanced
Simple feature collection with 431 features and 29 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 280371.1 ymin: 5235856 xmax: 921292.4 ymax: 6101487
Projected CRS: ETRS89 / UTM zone 32N
# A tibble: 431 × 30
     ADE    GF   BSG ARS   AGS   SDV_ARS     GEN   BEZ     IBZ BEM   NBD   SN_L 
   <int> <int> <int> <chr> <chr> <chr>       <chr> <chr> <int> <chr> <chr> <chr>
 1     4     4     1 01001 01001 0100100000… Flen… Krei…    40 --    ja    01   
 2     4     4     1 01002 01002 0100200000… Kiel  Krei…    40 --    ja    01   
 3     4     4     1 01003 01003 0100300000… Lübe… Krei…    40 --    ja    01   
 4     4     4     1 01004 01004 0100400000… Neum… Krei…    40 --    ja    01   
 5     4     4     1 01051 01051 0105100440… Dith… Kreis    42 --    ja    01   
 6     4     4     1 01053 01053 0105301001… Herz… Kreis    42 --    ja    01   
 7     4     4     1 01054 01054 0105400560… Nord… Kreis    42 --    ja    01   
 8     4     4     1 01055 01055 0105500120… Osth… Kreis    42 --    ja    01   
 9     4     4     1 01056 01056 0105600390… Pinn… Kreis    42 --    ja    01   
10     4     4     1 01057 01057 0105700570… Plön  Kreis    42 --    ja    01   
# ℹ 421 more rows
# ℹ 18 more variables: SN_R <chr>, SN_K <chr>, SN_V1 <chr>, SN_V2 <chr>,
#   SN_G <chr>, FK_S3 <chr>, NUTS <chr>, ARS_0 <chr>, AGS_0 <chr>, WSK <date>,
#   DEBKG_ID <chr>, geometry <MULTIPOLYGON [m]>, car_density <dbl>,
#   ecar_share <chr>, publictransport_meandist <dbl>, population <dbl>,
#   green_voteshare <dbl>, afd_voteshare <dbl>
echarging_sf
Simple feature collection with 60549 features and 5 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 5.243745 ymin: 47.2844 xmax: 15.54381 ymax: 55.50014
Geodetic CRS:  WGS 84
# A tibble: 60,549 × 6
   operator     federal_state power_kw type  num_plugs            geometry
 * <chr>        <chr>            <dbl> <chr>     <dbl>         <POINT [°]>
 1 deer GmbH    Baden-Württe…       44 Norm…         2  (9.72289 48.32789)
 2 EnBW mobili… Baden-Württe…       93 Schn…         2  (9.87484 48.57853)
 3 SWU Energie… Baden-Württe…       44 Norm…         2  (10.1934 48.52898)
 4 SWU Energie… Baden-Württe…       44 Norm…         2 (10.08268 48.55354)
 5 SWU Energie… Baden-Württe…       22 Norm…         1 (10.07698 48.17996)
 6 EnBW mobili… Baden-Württe…       30 Norm…         2 (9.980588 48.48039)
 7 SWU Energie… Baden-Württe…       22 Norm…         2 (9.985855 48.48316)
 8 Albwerk Gmb… Baden-Württe…       22 Norm…         2 (9.760458 48.46448)
 9 EnBW ODR AG  Baden-Württe…       22 Norm…         1 (10.02912 48.49882)
10 Autohaus Bu… Baden-Württe…       22 Norm…         2 (9.774676 48.40366)
# ℹ 60,539 more rows

We conduct a spatial join!

Adding district information to the point layer

echarging_sf_joined <-
  echarging_sf |> 
  sf::st_transform(
    crs = sf::st_crs(german_districts)
  ) |> 
  sf::st_join(
    german_districts |> 
      dplyr::select(AGS)
  )

plot(echarging_sf_joined["AGS"])

Spatial join even easier

Say, we want to count the number of charging stations in each German district.

# adjust crs first
echarging_sf <-
  sf::st_transform(echarging_sf, crs = sf::st_crs(german_districts))

# count all chargers in district
german_districts_enhanced$charger_in_districts <-
  lengths(sf::st_intersects(german_districts_enhanced, echarging_sf))

german_districts_enhanced$charger_in_districts
  [1]   67  208   95   94  156  106  259  215  180   72  239  185  169   56  162
 [16] 1223  216  121  452   92  100   48  101   69   36  199 1032  122   54  118
 [31]   33  111  117   77   63  114   26   98   78   96  161  136   48  130   33
 [46]  107   75  171   28   51  116  111  286   80  108   98   71  199  131   33
 [61]   28  342   46  771  136  348   57  144   27  120   65   79  119  204  366
 [76]  292  204  279  240  659   55  543  178  213  103  143  138  138  347   61
 [91]  105  167  361  155  240  219  140  157  220  125   95  173  160  279  297
[106]  563   88  102   59  251  201  218   69  119  197  294  107  474   37  153
[121]  173  191  535  145  329  145   67  235   98  157  146  145  125  108   86
[136]  174  252   57  217   88   98   43   88   74   46   81   23   52  123   91
[151]   84   40  115   66   94   60   55   79   27   71   61   79  120   48   21
[166]   45   36   45   56   71   23   40   41   22   57   88  143   35 1304  768
[181]  457  241  391  353  217  433  156  185  105   84  232   55  202  293  209
[196]  132  199  123  306   67  173  128   96  233  194   83  383   92  129  109
[211]  175  144  131  229  244  114  231  434   96  246  337   72  782 1634   28
[226]  120   91   93  102  105  115  127  173  150  111   89  144   60  656   88
[241]  113  271  107  153  117   57   80   44  121  125  104  177  200   72   69
[256]   82  396   12  423   32   58   86  114   33  202  119   38   75   46   37
[271]   17   84   51   49   88   54   42   50   36   55   26   89  131  392   14
[286]  181   98   81  126   85  128   64   28   67  124  103   51   55   43  136
[301]   35  106  108   95  165   40   69   46   91  157   51   79  136   85  102
[316]  123  110  197  149   59   44   92   80   42 2534   33   38   20   81   96
[331]   88   32   70   62  114   42  608   73  136   26   45  113   73   81   40
[346]  115  115  192   73  108   80  183  159  123  114  234  406  147  101   95
[361]  110  371  129   82   40   99   93   30   57   89   91   91   26   72   82
[376]   82   64   30  107   56   53   19   56   27   63   43   32   44   28   69
[391]   54   32   43   93   42   33   50   63   66   29   41    0    0    0    0
[406]    0    0    0    0    0    0    0    0    0    0    0    0    0    0    0
[421]    0    0    0    0    0    0    0    0    1    0    0

Charger count per district

tm_shape(german_districts_enhanced) +
  tm_polygons(
    fill = "charger_in_districts",
    fill.legend = 
      tm_legend(
        title = "Charger Count (Quantiles)"
      ),
    fill.scale = 
      tm_scale_intervals(
        style = "quantile", 
        values = "viridis"
      ),
    col = "lightgrey"
  )

3. Subsetting

Subsetting vector data

One might be interested in only one specific area of Germany, like Cologne. To subset a sf object, you can often use your usual data wrangling workflow. In this case, I know the AGS ID, and that is the only row I want to keep.

# subsetting
cologne <-
  german_districts_enhanced |> 
  dplyr::filter(AGS == "05315") |>  
  dplyr::select(AGS) 

plot(cologne)

Using sf for subsetting

If you have no information about ids but only about the geolocation, you can use sf::st_touches() (or sf::st_within(), sf::st_intersect(), sf::st_crosses(), …) to identify, for example, all districts which share a border with Cologne.

cologne_surrounding <-
  german_districts_enhanced |> 
  dplyr::select(AGS) |>  
  # length of mutual border > 0
  dplyr::filter(
    lengths(
      sf::st_touches(german_districts_enhanced, cologne)
    ) > 0
  ) 

plot(cologne_surrounding)

Cropping data

We can use our raster data to subset our e-charging data. But first, we have to adjust the CRS again (use sf::st_crs(echarging_sf) to look up EPSG code):

inhabitants_cologne <- terra::rast("./data/inhabitants_cologne.tif")

inhabitants_cologne <- terra::project(inhabitants_cologne, "EPSG:25832")

Now we can crop the data:

echarging_sf_cologne <-
  echarging_sf |> 
  sf::st_crop(inhabitants_cologne)

E-charging stations in Cologne

tm_shape(echarging_sf_cologne) +
  tm_dots()

‘Subsetting’ Raster Layers

As you know, we can subset vector data by simply filtering for specific attribute values. For example, to subset Cologne’s districts only by the one of Deutz, we can use the dplyr for sf data:

deutz <-
  sf::st_read(
    "./data/cologne.shp",
    quiet = TRUE
  ) |> 
  dplyr::filter(NAME == "Deutz") |> 
  sf::st_transform(25832)

tm_shape(deutz) +
  tm_borders()

Cropping raster data

Cropping is then a method of cutting out a specific slice of a raster layer based on an input dataset or geospatial extent, such as a bounding box.

cropped_inhabitants_cologne <-
  terra::crop(
    inhabitants_cologne, 
    deutz
  )

tm_shape(cropped_inhabitants_cologne) +
  tm_raster()

Masking

Masking is similar to cropping, yet values outside the extent are set to missing values (NA).

masked_inhabitants_cologne <-
  raster::mask(
    inhabitants_cologne, 
    deutz
  )

tm_shape(masked_inhabitants_cologne) +
  tm_raster()

Combining Cropping and Masking

cropped_masked_inhabitants_cologne <-
  terra::crop(
    inhabitants_cologne, 
    deutz
  ) |> 
  terra::mask(deutz)

tm_shape(cropped_masked_inhabitants_cologne) +
  tm_raster()

4. Extraction & Aggregation

Changes in terminology

In case, we only want to add one attribute from a vector dataset Y to another vector dataset X we can conduct a spatial join using sf::st_join() as before. There is nothing new to tell. In the raster data world, these kind of operations are called raster extractions.

Extracting information from raster data

Raster data are helpful when we aim to

  • Apply calculations that are the same for all geometries in the dataset
  • Extract information from raster fast and efficient

Do you remember these data operations?

immigrants_cologne <-  terra::rast("./data/immigrants_cologne.tif")

inhabitants_cologne <- terra::rast("./data/inhabitants_cologne.tif")

immigrant_rate <-
  immigrants_cologne * 100 / 
  inhabitants_cologne

Raster extraction

To extract the raster values at a specific point by location, we use the following:

terra::extract(immigrant_rate, echarging_sf_cologne, ID = FALSE)
    immigrants_cologne
1                   NA
2                   NA
3                   NA
4                   NA
5                   NA
6                   NA
7                   NA
8                   NA
9                   NA
10                  NA
11           25.274725
12           25.274725
13                  NA
14                  NA
15           19.444444
16                  NA
17                  NA
18                  NA
19                  NA
20                  NA
21                  NA
22                  NA
23                  NA
24           10.000000
25                  NA
26                  NA
27                  NA
28                  NA
29                  NA
30            5.389222
31                  NA
32                  NA
33           12.790698
34            8.737864
35            4.705882
36                  NA
37                  NA
38                  NA
39                  NA
40                  NA
41                  NA
42                  NA
43                  NA
44                  NA
45                  NA
46                  NA
47           12.500000
48                  NA
49                  NA
50           29.729730
51           50.000000
52                  NA
53            8.333333
54            8.333333
55           12.142857
56                  NA
57                  NA
58                  NA
59            5.633803
60           17.647059
61            8.000000
62                  NA
63                  NA
64           10.526316
65           16.883117
66           21.739130
67           50.000000
68                  NA
69                  NA
70                  NA
71                  NA
72                  NA
73                  NA
74                  NA
75                  NA
76                  NA
77                  NA
78                  NA
79                  NA
80                  NA
81                  NA
82                  NA
83                  NA
84                  NA
85                  NA
86                  NA
87                  NA
88            9.677419
89            3.809524
90                  NA
91                  NA
92                  NA
93                  NA
94                  NA
95                  NA
96                  NA
97                  NA
98                  NA
99                  NA
100                 NA
101          15.000000
102           7.894737
103                 NA
104          18.000000
105                 NA
106                 NA
107          11.904762
108          16.666667
109          13.636364
110                 NA
111                 NA
112                 NA
113          17.142857
114           6.818182
115                 NA
116                 NA
117                 NA
118          17.142857
119          10.526316
120                 NA
121                 NA
122          33.043478
123                 NA
124          13.888889
125                 NA
126                 NA
127          30.136986
128                 NA
129                 NA
130          26.086957
131           7.407407
132           7.500000
133                 NA
134          14.285714
135                 NA
136          75.000000
137                 NA
138          15.789474
139           7.692308
140          10.828025
141                 NA
142                 NA
143                 NA
144                 NA
145                 NA
146                 NA
147                 NA
148                 NA
149                 NA
150                 NA
151                 NA
152                 NA
153                 NA
154                 NA
155                 NA
156                 NA
157                 NA
158                 NA
159                 NA
160                 NA
161           9.090909
162                 NA
163                 NA
164                 NA
165                 NA
166                 NA
167                 NA
168                 NA
169                 NA
170          29.729730
171                 NA
172                 NA
173                 NA
174                 NA
175                 NA
176          17.857143
177          18.309859
178          18.309859
179          18.309859
180          18.309859
181          18.309859
182          18.309859
183          13.043478
184          13.043478
185          13.043478
186          13.043478
187          13.043478
188          13.043478
189          13.043478
190          13.043478
191          13.043478
192           6.629834
193          23.076923
194          18.181818
195          23.076923
196          23.076923
197          17.582418
198           7.951807
199                 NA
200          16.858238
201          18.614719
202          16.814159
203          17.391304
204          21.904762
205          18.652850
206           9.375000
207           9.375000
208                 NA
209                 NA
210                 NA
211                 NA
212                 NA
213                 NA
214                 NA
215                 NA
216                 NA
217                 NA
218                 NA
219                 NA
220                 NA
221                 NA
222                 NA
223                 NA
224                 NA
225                 NA
226                 NA
227                 NA
228                 NA
229                 NA
230                 NA
231                 NA
232                 NA
233                 NA
234                 NA
235                 NA
236          10.582011
237          18.072289
238           8.679245
239          27.083333
240          11.976048
241          18.500000
242          13.281250
243          13.809524
244           7.594937
245           6.611570
246          12.195122
247          19.780220
248                 NA
249          12.871287
250          15.384615
251          16.438356
252          22.857143
253           4.347826
254                 NA
255          11.940299
256                 NA
257          51.063830
258          20.720721
259          10.924370
260                 NA
261                 NA
262                 NA
263                 NA
264                 NA
265                 NA
266                 NA
267                 NA
268                 NA
269                 NA
270                 NA
271                 NA
272                 NA
273                 NA
274                 NA
275                 NA
276                 NA
277                 NA
278                 NA
279                 NA
280                 NA
281                 NA
282                 NA
283                 NA
284                 NA
285                 NA
286                 NA
287          20.600858
288          15.636364
289                 NA
290                 NA
291          25.373134
292          11.956522
293                 NA
294                 NA
295                 NA
296          25.773196
297          14.942529
298          14.534884
299           9.090909
300          30.400000
301          26.737968
302          15.853659
303                 NA
304          20.000000
305          29.166667
306                 NA
307          11.428571
308          27.192982
309                 NA
310                 NA
311                 NA
312          19.473684
313           3.846154
314           9.090909
315          48.275862
316          14.110429
317                 NA
318          14.210526
319                 NA
320          26.315789
321          18.750000
322           6.000000
323                 NA
324           6.382979
325           6.382979
326           8.602151
327           8.602151
328                 NA
329                 NA
330          33.834586
331                 NA
332          26.785714
333           5.660377
334          29.411765
335                 NA
336                 NA
337          22.142857
338                 NA
339                 NA
340          26.885880
341                 NA
342          30.303030
343                 NA
344                 NA
345                 NA
346          13.043478
347           8.333333
348           8.333333
349          17.391304
350           7.878788
351                 NA
352                 NA
353                 NA
354           5.714286
355          28.828829
356                 NA
357                 NA
358                 NA
359          12.903226
360          12.000000
361                 NA
362                 NA
363                 NA
364                 NA
365                 NA
366                 NA
367                 NA
368                 NA
369                 NA
370                 NA
371                 NA
372                 NA
373                 NA
374          26.315789
375          10.900474
376          28.571429
377          27.513228
378          27.513228
379           8.677686
380          24.686192
381                 NA
382                 NA
383                 NA
384                 NA
385                 NA
386                 NA
387                 NA
388                 NA
389                 NA
390                 NA
391                 NA
392                 NA
393                 NA
394                 NA
395                 NA
396                 NA
397                 NA
398                 NA
399                 NA
400                 NA
401                 NA
402                 NA
403                 NA
404                 NA
405          20.574163
406           2.684564
407          19.469027
408          10.077519
409                 NA
410                 NA
411          12.418301
412          23.529412
413                 NA
414          15.000000
415           9.000000
416           8.943089
417          25.581395
418          25.581395
419          32.142857
420                 NA
421                 NA
422                 NA
423                 NA
424                 NA
425                 NA
426                 NA
427                 NA
428                 NA
429          16.603774
430                 NA
431                 NA
432                 NA
433                 NA
434                 NA
435                 NA
436                 NA
437                 NA
438                 NA
439                 NA
440                 NA
441                 NA
442                 NA
443                 NA
444                 NA
445                 NA
446                 NA
447                 NA
448                 NA
449                 NA
450                 NA
451                 NA
452                 NA
453                 NA
454                 NA
455                 NA
456                 NA
457                 NA
458                 NA
459                 NA
460                 NA
461                 NA
462                 NA
463                 NA
464                 NA
465                 NA
466                 NA
467                 NA
468                 NA
469                 NA
470                 NA
471                 NA
472                 NA
473                 NA
474                 NA
475                 NA
476                 NA
477                 NA
478                 NA
479                 NA
480                 NA
481                 NA
482                 NA
483                 NA
484                 NA
485                 NA
486                 NA
487                 NA
488                 NA
489                 NA
490                 NA
491                 NA
492                 NA
493                 NA
494                 NA
495                 NA
496                 NA
497                 NA
498                 NA
499                 NA
500                 NA
501                 NA
502                 NA
503                 NA
504                 NA
505                 NA
506                 NA
507                 NA
508                 NA
509                 NA
510                 NA
511                 NA
512                 NA
513                 NA
514                 NA
515                 NA
516                 NA
517                 NA
518                 NA
519                 NA
520                 NA
521                 NA
522                 NA
523                 NA
524                 NA
525                 NA
526                 NA
527                 NA
528                 NA
529                 NA
530                 NA
531                 NA
532                 NA
533                 NA
534           8.823529
535                 NA
536                 NA
537                 NA
538                 NA
539                 NA
540                 NA
541                 NA
542                 NA
543                 NA
544                 NA
545                 NA
546                 NA
547           6.521739
548                 NA
549          16.279070
550          10.169492
551                 NA
552                 NA
553                 NA
554                 NA
555          57.142857
556                 NA
557                 NA
558                 NA
559                 NA
560                 NA
561          18.320611
562                 NA
563           9.883721
564                 NA
565           6.382979
566                 NA
567                 NA
568                 NA
569                 NA
570                 NA
571                 NA
572                 NA
573                 NA
574                 NA
575                 NA
576                 NA
577           8.695652
578          60.000000
579          15.079365
580                 NA
581           7.272727
582           6.818182
583           4.761905
584           3.333333
585           6.481481
586           6.285714
587          21.666667
588          10.416667
589           7.377049
590          10.526316
591                 NA
592                 NA
593                 NA
594          21.276596
595           5.042017
596          10.650888
597          10.169492
598                 NA
599           3.703704
600          43.767313
601                 NA
602                 NA
603                 NA
604                 NA
605                 NA
606                 NA
607                 NA
608           8.000000
609          20.833333
610           5.555556
611          20.714286
612          12.121212
613          11.688312
614          23.809524
615          13.043478
616          23.004695
617                 NA
618                 NA
619                 NA
620          24.324324
621           4.278075
622          15.131579
623          12.658228
624          24.742268
625          45.833333
626                 NA
627           9.779180
628          11.666667
629                 NA
630                 NA
631          21.428571
632          21.428571
633          21.428571
634          21.428571
635           6.097561
636          26.595745
637                 NA
638          12.121212
639                 NA
640                 NA
641                 NA
642          26.114650
643          33.014354
644          54.237288
645          54.237288
646          19.444444
647          23.809524
648          50.000000
649          35.256410
650                 NA
651          41.007194
652          11.267606
653          23.809524
654                 NA
655          16.666667
656                 NA
657          21.717172
658          20.670391
659          21.243523
660                 NA
661          13.402062
662          28.125000
663           6.666667
664           6.666667
665                 NA
666          20.661157
667           4.166667
668          23.076923
669          20.353982
670          16.000000
671          10.000000
672          11.607143
673          30.909091
674                 NA
675                 NA
676                 NA
677                 NA
678                 NA
679                 NA
680                 NA
681                 NA
682                 NA
683                 NA
684                 NA
685                 NA
686                 NA
687                 NA
688                 NA
689                 NA
690                 NA
691          30.769231
692          20.430108
693          10.526316
694          46.590909
695          15.463918
696                 NA
697          18.390805
698           5.217391
699                 NA
700           5.555556
701           2.678571
702          37.404580
703                 NA
704          26.395939
705                 NA
706           5.555556
707           6.818182
708           4.347826
709           4.347826
710          33.043478
711          43.859649
712          33.039648
713           2.803738
714           7.500000
715          11.904762
716                 NA
717                 NA
718          32.467532
719           7.633588
720          13.861386
721          14.285714
722                 NA
723          22.222222
724          12.643678
725                 NA
726                 NA
727                 NA
728           4.761905
729          53.846154
730                 NA
731                 NA
732                 NA
733          42.424242
734          11.111111
735           4.477612
736           6.250000
737                 NA
738                 NA
739                 NA
740          42.857143
741                 NA
742                 NA
743                 NA
744                 NA
745                 NA
746                 NA
747                 NA
748                 NA
749                 NA
750                 NA
751          16.000000
752                 NA
753                 NA
754                 NA
755                 NA
756                 NA
757                 NA
758                 NA
759                 NA
760                 NA
761                 NA
762                 NA
763                 NA
764                 NA
765                 NA
766                 NA
767                 NA
768                 NA
769          24.844720
770                 NA
771          20.500000
772                 NA
773                 NA
774          31.034483
775                 NA
776                 NA
777                 NA
778                 NA
779                 NA
780                 NA
781                 NA
782                 NA
783                 NA
784                 NA
785                 NA
786                 NA
787                 NA
788                 NA
789                 NA
790          18.750000
791                 NA
792                 NA
793          18.750000
794          18.750000
795                 NA
796                 NA
797          60.000000
798                 NA
799                 NA
800                 NA
801                 NA
802                 NA
803                 NA
804                 NA
805                 NA
806                 NA
807                 NA
808                 NA
809                 NA
810                 NA
811                 NA
812                 NA
813                 NA
814                 NA
815                 NA
816          25.925926
817                 NA
818           7.563025
819                 NA
820          61.111111
821                 NA
822                 NA
823           8.000000
824                 NA
825          43.023256
826                 NA
827                 NA
828          20.289855
829          21.739130
830                 NA
831                 NA
832                 NA
833          62.500000
834                 NA
835                 NA
836                 NA
837                 NA
838                 NA
839          28.205128
840                 NA
841          17.307692
842          12.837838
843          15.625000
844          27.741935
845          15.596330
846          42.000000
847          13.043478
848          48.648649
849                 NA
850          20.833333
851          20.833333
852           7.936508
853                 NA
854                 NA
855           4.838710
856           4.838710
857                 NA
858                 NA
859                 NA
860                 NA
861                 NA
862                 NA
863                 NA
864                 NA
865          29.032258
866                 NA
867          10.000000
868                 NA
869                 NA
870          27.027027
871                 NA
872                 NA
873           6.976744
874                 NA
875                 NA
876           4.687500

Add results to existing dataset

This information can be added to an existing dataset (our points in this example):

echarging_sf_cologne$immigrant_rate_value <-
  terra::extract(immigrant_rate, echarging_sf_cologne, ID = FALSE) |> 
  unlist()

echarging_sf_cologne
Simple feature collection with 876 features and 6 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 343957.8 ymin: 5632513 xmax: 370664.4 ymax: 5661644
Projected CRS: ETRS89 / UTM zone 32N
# A tibble: 876 × 7
   operator     federal_state power_kw type  num_plugs           geometry
 * <chr>        <chr>            <dbl> <chr>     <dbl>        <POINT [m]>
 1 Cut! Energy… Nordrhein-We…     44   Norm…         2 (358799.3 5659476)
 2 Kaufland     Nordrhein-We…     93   Schn…         2 (352129.2 5661644)
 3 Kaufland     Nordrhein-We…     72   Schn…         3 (352135.4 5661628)
 4 MEGA Kunden… Nordrhein-We…     44   Norm…         2 (352102.1 5661543)
 5 MEGA Kunden… Nordrhein-We…     44   Norm…         2 (352102.1 5661543)
 6 MEGA Kunden… Nordrhein-We…     44   Norm…         2 (352102.1 5661543)
 7 MEGA Kunden… Nordrhein-We…     44   Norm…         2 (352102.1 5661543)
 8 MEGA Kunden… Nordrhein-We…     44   Norm…         2 (352102.1 5661543)
 9 Auto Thomas… Nordrhein-We…     39.6 Norm…         2 (352737.2 5633750)
10 EWE Go GmbH  Nordrhein-We…    160   Schn…         2   (353129 5634718)
# ℹ 866 more rows
# ℹ 1 more variable: immigrant_rate_value <dbl>

More elaborated: spatial buffers

Sometimes, extracting information 1:1 is not enough

  • Too narrow
  • Missing information about the surroundings of a point
tm_shape(immigrant_rate) +
  tm_raster() +
  tm_shape(
    sf::st_buffer(echarging_sf_cologne, 500) 
  ) +
  tm_dots(size = .1) +
  tm_borders()

Buffer extraction

We can use spatial buffers of different sizes to extract information on surroundings:

terra::extract(
  immigrant_rate, 
  sf::st_buffer(echarging_sf_cologne, 500),
  fun = mean,
  na.rm = TRUE,
  ID = FALSE,
  raw = TRUE
)
    immigrants_cologne
1                  NaN
2            30.689268
3            31.024269
4            30.109133
5            30.109133
6            30.109133
7            30.109133
8            30.109133
9            18.839996
10                 NaN
11           18.251163
12           18.251163
13           12.480011
14           12.480011
15           18.424806
16                 NaN
17                 NaN
18                 NaN
19           19.178225
20                 NaN
21                 NaN
22                 NaN
23                 NaN
24           10.514073
25           62.500000
26           62.500000
27           62.500000
28          100.000000
29           13.060481
30           13.757629
31           23.035654
32           23.035654
33           12.874063
34           15.714064
35           17.759908
36           15.964036
37           17.274671
38           22.558999
39                 NaN
40                 NaN
41                 NaN
42                 NaN
43                 NaN
44                 NaN
45           25.203057
46           21.455188
47           21.541857
48           20.330240
49           20.330240
50           17.108656
51           17.459820
52           43.825944
53           13.727532
54           13.727532
55           13.341282
56           16.776713
57           21.047340
58            9.702687
59           11.442805
60           13.724539
61           13.322034
62           20.663649
63           11.407075
64           10.034055
65           15.578801
66           16.972791
67           19.161427
68           26.867816
69           26.867816
70           26.867816
71           26.867816
72           26.867816
73            9.364526
74           13.405072
75           13.620435
76           27.272727
77           18.133372
78           14.677015
79           14.677015
80                 NaN
81           15.000000
82          100.000000
83          100.000000
84          100.000000
85           14.779270
86           14.779270
87            9.876248
88            9.609889
89            9.961877
90           22.542760
91           23.789690
92           24.123839
93           20.684094
94           23.367645
95           14.572866
96           12.336835
97           14.275269
98           12.928160
99           12.928160
100           9.126074
101           9.354862
102          11.965110
103          12.511551
104           9.490522
105           9.071518
106           9.579474
107          10.778239
108           8.275924
109          13.897536
110          16.070802
111          11.480752
112          12.928155
113          32.536811
114          19.039160
115          24.478408
116          31.939751
117          31.939751
118          32.536811
119          30.280234
120          24.863629
121          24.863629
122          25.069421
123          22.289562
124          10.165677
125          14.443126
126          22.156499
127          22.014963
128          27.883306
129          14.279433
130          21.438430
131          11.816043
132          13.755265
133          21.246422
134          23.110947
135          22.167128
136          28.072219
137          17.239800
138          18.747400
139          17.423669
140          19.247369
141           7.436713
142           7.436713
143          13.106685
144          13.106685
145                NaN
146          15.188014
147          19.565217
148          15.188014
149          16.348265
150          15.188014
151          16.348265
152          15.188014
153          16.348265
154          16.348265
155          20.286358
156          21.543604
157          21.543604
158          14.376509
159           8.092493
160           7.589144
161          15.637289
162          13.999447
163          22.226983
164           7.170879
165                NaN
166          34.848485
167          34.848485
168          34.848485
169          34.848485
170          16.690038
171          19.077622
172          14.636317
173          14.636317
174          14.636317
175          14.636317
176          17.517823
177          18.121414
178          18.121414
179          18.121414
180          18.121414
181          18.121414
182          18.121414
183          18.057020
184          18.057020
185          18.242471
186          18.156316
187          18.242471
188          18.242471
189          18.242471
190          18.242471
191          18.242471
192          20.322962
193          24.603270
194          25.726868
195          24.603270
196          24.558092
197          26.255361
198          10.362776
199          10.022889
200          13.157760
201          23.463178
202          12.282973
203          20.313225
204          24.837026
205          18.532695
206          17.517177
207          17.517177
208          17.555169
209          17.555169
210          17.555169
211          17.555169
212          17.555169
213          17.555169
214          17.555169
215          17.555169
216          17.555169
217          17.555169
218          17.555169
219          17.555169
220          17.555169
221          17.555169
222          17.555169
223          17.555169
224          17.555169
225          17.555169
226          17.555169
227          17.555169
228          17.555169
229          17.555169
230          17.555169
231          17.555169
232          17.555169
233          17.555169
234          17.555169
235          15.317597
236          16.217751
237          13.972405
238          15.496273
239          16.670677
240          15.186751
241          15.217877
242          14.818024
243          14.719214
244          14.098589
245          21.577599
246          17.093786
247          16.779919
248          17.537970
249          15.139761
250          14.871517
251          26.271610
252          14.072582
253          14.129258
254          12.594627
255          16.989394
256          15.727774
257          19.089677
258          16.285333
259          15.900499
260          35.290712
261          35.290712
262          35.290712
263          35.290712
264          24.772475
265          24.772475
266          32.329089
267          15.789474
268          15.789474
269          18.357558
270          18.357558
271          34.010199
272          34.010199
273          34.010199
274          34.010199
275          34.010199
276          34.010199
277          34.010199
278          34.010199
279          34.010199
280          34.010199
281          12.840996
282          12.840996
283          12.840996
284          12.840996
285          12.840996
286          12.840996
287          12.239287
288          12.760854
289          14.702236
290          12.626195
291          34.215263
292          17.230151
293          14.940545
294          14.940545
295          17.363913
296          16.340990
297          16.063010
298          16.010140
299          15.288145
300          19.480073
301          16.506619
302          17.654480
303          17.456971
304          15.032626
305          20.492259
306          16.360149
307          15.624051
308          20.481348
309          23.464551
310          22.438974
311          22.438974
312          21.276037
313          17.000543
314          17.366188
315          23.981567
316          25.120230
317          84.848485
318          13.646955
319          11.321193
320          21.979540
321          15.152983
322           9.706929
323          56.708683
324           9.971849
325           9.971849
326          19.627202
327          19.627202
328          32.470658
329          31.198691
330          23.039739
331          47.236467
332          15.249318
333          20.520565
334          28.675118
335          20.509425
336          47.236467
337          23.357945
338          39.214676
339          18.936742
340          21.286377
341           8.266507
342          18.764553
343          14.757196
344          13.064173
345          13.064173
346          13.692316
347          15.401425
348          15.401425
349          10.875724
350          15.824836
351           8.451705
352           8.451705
353          11.524785
354          10.039274
355          17.552924
356          65.789474
357          65.789474
358          13.254952
359          21.114722
360          17.029447
361          30.673487
362          27.013731
363          10.861871
364          10.861871
365          10.861871
366          10.861871
367          10.861871
368          10.861871
369          10.861871
370          10.861871
371          10.861871
372          10.861871
373          10.861871
374          13.708125
375          17.812532
376          18.608512
377          15.550220
378          15.550220
379          18.913155
380          16.173418
381          31.149886
382          31.149886
383          31.149886
384          24.305103
385          24.305103
386          28.726909
387          23.200644
388          23.200644
389          24.186564
390          24.186564
391          24.186564
392          24.186564
393          24.186564
394          24.186564
395          24.186564
396          24.186564
397          24.186564
398          31.212990
399          31.212990
400          23.820456
401          23.820456
402          23.820456
403          23.820456
404          23.820456
405          24.060261
406           8.534757
407          13.704405
408          10.700428
409          14.099080
410          25.664434
411          18.785299
412          17.436260
413          20.345449
414          20.844770
415          15.386900
416          17.168522
417          14.320328
418          14.320328
419          17.203823
420                NaN
421                NaN
422                NaN
423                NaN
424                NaN
425                NaN
426                NaN
427          23.233950
428          23.233950
429          20.182639
430          23.964897
431          23.964897
432         100.000000
433         100.000000
434         100.000000
435         100.000000
436         100.000000
437         100.000000
438         100.000000
439         100.000000
440         100.000000
441         100.000000
442         100.000000
443         100.000000
444         100.000000
445         100.000000
446         100.000000
447         100.000000
448         100.000000
449         100.000000
450         100.000000
451         100.000000
452         100.000000
453         100.000000
454         100.000000
455         100.000000
456         100.000000
457         100.000000
458         100.000000
459         100.000000
460         100.000000
461         100.000000
462         100.000000
463         100.000000
464         100.000000
465         100.000000
466         100.000000
467         100.000000
468         100.000000
469         100.000000
470         100.000000
471         100.000000
472         100.000000
473         100.000000
474         100.000000
475         100.000000
476         100.000000
477         100.000000
478         100.000000
479         100.000000
480         100.000000
481         100.000000
482         100.000000
483         100.000000
484         100.000000
485         100.000000
486         100.000000
487         100.000000
488         100.000000
489         100.000000
490         100.000000
491         100.000000
492         100.000000
493         100.000000
494                NaN
495          11.113062
496                NaN
497                NaN
498                NaN
499                NaN
500                NaN
501                NaN
502                NaN
503                NaN
504                NaN
505                NaN
506                NaN
507                NaN
508                NaN
509                NaN
510                NaN
511                NaN
512                NaN
513                NaN
514                NaN
515                NaN
516                NaN
517                NaN
518                NaN
519                NaN
520                NaN
521                NaN
522                NaN
523                NaN
524                NaN
525                NaN
526                NaN
527                NaN
528                NaN
529                NaN
530                NaN
531                NaN
532                NaN
533          13.727655
534          13.354236
535          60.864662
536          44.772098
537          39.723602
538          39.723602
539          98.076923
540          98.076923
541          62.105263
542          62.105263
543          62.105263
544          62.105263
545          62.105263
546          62.105263
547           9.353183
548          13.993590
549          11.253390
550           8.894746
551          62.105263
552          13.159557
553          13.159557
554          13.159557
555          11.846289
556          26.818880
557          19.578483
558          11.859017
559          11.859017
560           9.242620
561          10.158558
562          11.958980
563          10.790331
564          10.525738
565           9.741421
566          10.649967
567          19.627231
568          19.627231
569          53.846154
570           8.120422
571          22.043419
572          22.043419
573          22.043419
574          22.043419
575          13.366921
576          12.777159
577          16.252081
578          20.549149
579          23.618988
580          19.867754
581          10.228215
582           6.255870
583           6.256628
584           7.117795
585          11.789633
586           6.089852
587          10.245778
588           9.135267
589           8.969029
590           7.427639
591          10.294174
592           7.026119
593          14.949885
594          13.872839
595           8.432926
596           7.369677
597           7.989067
598           7.181416
599           7.836296
600          13.415633
601          12.715183
602          12.715183
603          12.715183
604          12.715183
605          19.451122
606          13.480584
607          16.344883
608          18.629271
609          12.235893
610          11.406428
611          18.146996
612          13.100882
613          12.259738
614          16.268679
615          13.454494
616          15.301908
617          13.914263
618          15.415713
619          17.633649
620          13.269333
621          14.579759
622          14.239716
623          14.361029
624          20.448234
625          14.823154
626          65.000000
627          12.847498
628          10.303374
629                NaN
630                NaN
631          18.129126
632          18.129126
633          18.129126
634          18.129126
635          10.212428
636          10.661272
637           9.031977
638           8.800828
639          12.066649
640          12.787486
641          15.038232
642          14.871538
643          15.780669
644          40.995471
645          40.995471
646          21.265517
647          29.175205
648          31.191287
649          34.186923
650          28.196166
651          34.290973
652          23.547063
653          29.175205
654          22.028637
655          26.771218
656          28.546290
657          25.648111
658          22.185476
659          23.219860
660          31.865458
661          23.215181
662          23.062183
663          13.143755
664          13.143755
665          12.298931
666          12.690536
667          16.393625
668          19.010931
669          11.852314
670          18.370018
671          15.874996
672          10.040861
673          13.417263
674          14.145860
675          14.145860
676          14.145860
677          14.145860
678          14.145860
679          14.145860
680          14.145860
681          14.145860
682          14.145860
683          14.145860
684          12.454843
685          38.653298
686          38.653298
687          38.653298
688          38.653298
689          28.324210
690          34.714791
691          32.784513
692          26.079626
693          23.997821
694          33.036525
695          31.970639
696          40.297253
697          22.035236
698          12.657063
699          24.101173
700          12.884953
701          10.278415
702          30.000101
703          31.128699
704          29.330846
705          19.940401
706          12.884953
707           8.635181
708          13.289133
709          13.289133
710          22.956131
711          22.310157
712          30.173941
713          15.559987
714          17.300691
715          22.070435
716          17.885822
717          17.885822
718          10.605523
719          10.896273
720          25.140956
721          22.035657
722          22.484851
723          19.920496
724          11.275934
725          12.078559
726          12.005008
727          22.307356
728          10.849367
729          17.484287
730          13.566540
731          13.566540
732          17.121233
733          16.876633
734          15.377812
735          12.848985
736          13.699767
737           9.193501
738          16.669944
739          16.669944
740          32.142857
741                NaN
742                NaN
743                NaN
744                NaN
745                NaN
746                NaN
747                NaN
748                NaN
749                NaN
750                NaN
751          15.474401
752          36.363636
753          34.848485
754          17.612771
755          17.612771
756          17.612771
757          17.612771
758          56.565657
759          11.090226
760          11.090226
761          11.090226
762          11.090226
763          11.090226
764          11.090226
765          11.090226
766          11.090226
767          11.090226
768          11.090226
769          13.015290
770          56.565657
771          12.452092
772          31.740811
773          56.565657
774          31.920865
775          36.363636
776          36.363636
777          36.363636
778          36.363636
779          36.363636
780          36.363636
781          36.363636
782          36.363636
783         100.000000
784         100.000000
785          35.266683
786          34.247404
787          34.247404
788          34.247404
789          34.247404
790          34.260406
791          35.266683
792          35.266683
793          36.004230
794          37.731826
795          35.266683
796          34.214616
797          28.169265
798          31.833447
799          31.833447
800          32.660603
801          29.076582
802          30.932216
803          31.786471
804          30.932216
805          31.156846
806          31.156846
807          10.135115
808          10.135115
809          10.135115
810          10.135115
811          10.135115
812          10.135115
813          36.363636
814          36.363636
815          36.363636
816          21.184916
817          27.932543
818          15.584326
819          15.391452
820          14.481623
821          15.022787
822          26.179022
823          14.545675
824           9.627873
825          16.464436
826          16.197683
827          16.737836
828          14.678709
829          15.338412
830          15.622335
831          15.622335
832          15.622335
833          15.062795
834          17.432447
835          17.432447
836           9.232955
837          25.859851
838          25.859851
839          26.849386
840          16.705786
841          24.333115
842          24.680768
843          25.831182
844          22.823124
845          17.959382
846          26.821140
847          23.171492
848          30.924144
849          17.283298
850          15.362022
851          15.493390
852          10.046466
853          18.886424
854          13.333771
855           9.811799
856           9.888126
857          21.375721
858          22.540009
859          19.087613
860          19.087613
861          18.442048
862          18.442048
863          16.102033
864          14.746257
865          19.609447
866          15.998127
867          20.087950
868          20.583743
869          21.213259
870          22.658779
871          19.166508
872          15.654887
873          21.056426
874          19.530130
875          20.903429
876          15.252328
terra::extract(
  immigrant_rate, 
  sf::st_buffer(echarging_sf_cologne, 1000),
  fun = mean,
  na.rm = TRUE,
  ID = FALSE,
  raw = TRUE
)
    immigrants_cologne
1            13.529080
2            19.929001
3            19.607900
4            19.642150
5            19.642150
6            19.642150
7            19.642150
8            19.642150
9            18.302595
10           41.079433
11           14.025739
12           14.025739
13           10.286388
14           10.286388
15           15.970992
16           14.460422
17           14.460422
18           14.460422
19           16.177228
20           45.395123
21           45.395123
22           45.395123
23           45.395123
24           11.329216
25           32.173972
26           32.173972
27           32.173972
28           23.334825
29           15.152595
30           16.353176
31           17.166848
32           17.300654
33           12.400199
34           16.232756
35           14.641138
36           14.388859
37           15.886165
38           19.048553
39                 NaN
40                 NaN
41                 NaN
42                 NaN
43                 NaN
44                 NaN
45           15.960944
46           23.527919
47           24.629666
48           15.149836
49           15.149836
50           15.131865
51           15.362907
52           17.134482
53           13.662567
54           13.662567
55           15.898146
56           14.918616
57           16.938374
58           11.471347
59           11.729408
60           15.591820
61           14.140996
62           15.661559
63           10.980812
64           12.113624
65           15.113860
66           14.716975
67           15.507672
68           15.316302
69           15.316302
70           15.316302
71           15.316302
72           15.316302
73           10.119134
74           14.441562
75           15.003324
76           14.360505
77           13.135527
78           12.542610
79           12.542610
80           12.831605
81           15.193517
82           25.333156
83           25.333156
84           35.385552
85           16.410248
86           16.410248
87           12.385829
88           11.235300
89           12.267541
90           19.690377
91           19.690377
92           19.053614
93           19.690377
94           19.690377
95           16.565715
96           17.706569
97           15.269168
98           15.119310
99           15.119310
100           8.324908
101           8.856053
102          12.275663
103          13.382043
104           8.647714
105           8.402760
106           8.762963
107          13.040760
108           8.476846
109          13.449405
110          14.158421
111          13.732256
112          14.047107
113          21.967958
114          17.287955
115          18.973852
116          20.979583
117          20.979583
118          21.967958
119          23.137272
120          19.903420
121          19.903420
122          22.216930
123          20.029443
124          11.772214
125          10.861788
126          18.728154
127          19.509789
128          21.284422
129          14.398028
130          18.041181
131          12.779830
132          11.208671
133          17.597679
134          18.222336
135          21.284009
136          19.986349
137          20.048481
138          19.981407
139          20.187700
140          20.223886
141          13.656888
142          13.656888
143          14.886911
144          14.886911
145          15.896739
146          12.154482
147          13.073200
148          12.154482
149          13.812172
150          12.154482
151          13.812172
152          12.154482
153          13.812172
154          16.348265
155          21.464652
156          21.658955
157          21.658955
158          12.880951
159           8.281682
160           9.071311
161          17.325523
162          16.270742
163          18.408431
164          15.551314
165           7.110812
166          34.177627
167          34.177627
168          34.177627
169          34.177627
170          14.510839
171          16.507341
172          16.726791
173          16.726791
174          16.726791
175          16.726791
176          16.863387
177          17.447456
178          17.447456
179          17.447456
180          17.447456
181          17.447456
182          17.447456
183          18.483128
184          18.637448
185          18.637448
186          18.637448
187          18.598172
188          18.780315
189          18.780315
190          18.695913
191          18.598172
192          20.236247
193          20.588736
194          21.834017
195          20.588736
196          20.866374
197          22.286882
198          14.310459
199          13.675897
200          17.520898
201          21.838932
202          16.722790
203          18.882864
204          20.841537
205          21.386287
206          19.293722
207          19.293722
208          20.155043
209          20.155043
210          20.155043
211          20.155043
212          20.155043
213          20.155043
214          20.155043
215          20.155043
216          20.155043
217          20.155043
218          20.155043
219          20.155043
220          20.155043
221          20.155043
222          20.155043
223          20.155043
224          20.155043
225          20.155043
226          20.155043
227          20.155043
228          20.155043
229          20.155043
230          20.155043
231          20.155043
232          20.155043
233          20.155043
234          20.155043
235          16.781009
236          16.212651
237          17.233628
238          15.127043
239          13.691441
240          15.429402
241          14.248644
242          14.896604
243          14.579439
244          15.580156
245          19.953855
246          17.055723
247          16.733021
248          17.445029
249          16.281502
250          14.810606
251          21.399660
252          14.641620
253          15.855764
254          15.584692
255          15.972472
256          17.398875
257          15.733400
258          13.046227
259          16.194556
260          25.161658
261          25.161658
262          25.161658
263          25.161658
264          26.538657
265          26.538657
266          31.683341
267          17.385328
268          17.385328
269          18.695120
270          18.695120
271          32.924059
272          32.924059
273          32.924059
274          32.924059
275          32.924059
276          32.924059
277          32.924059
278          32.924059
279          32.924059
280          32.924059
281          14.312241
282          14.312241
283          14.312241
284          14.312241
285          14.312241
286          14.312241
287          14.199153
288          14.331740
289          15.520445
290          13.485136
291          29.449753
292          15.557177
293          16.496541
294          16.496541
295          17.023705
296          18.232128
297          16.643100
298          15.950695
299          17.992402
300          19.836864
301          16.689038
302          16.608271
303          18.820708
304          15.439799
305          17.396782
306          18.147558
307          18.241850
308          18.788584
309          21.868578
310          21.648835
311          21.648835
312          21.080038
313          19.481539
314          18.881549
315          22.136021
316          22.962252
317          64.976690
318          14.506221
319          11.587543
320          17.869369
321          17.246140
322          11.894196
323          31.774922
324          11.244498
325          11.244498
326          23.095500
327          23.095500
328          29.220682
329          30.978386
330          19.253464
331          34.473270
332          17.109129
333          22.288829
334          20.837800
335          19.629086
336          34.963890
337          19.482835
338          34.371841
339          18.822929
340          24.661054
341           9.579580
342          23.046852
343          12.464662
344          13.405247
345          13.405247
346          18.470715
347          13.117487
348          13.117487
349          12.162067
350          12.778851
351          12.729771
352          12.729771
353          17.017888
354          13.885893
355          15.552062
356          40.315291
357          44.249248
358          12.378056
359          26.264235
360          17.730731
361          21.695002
362          17.951643
363          17.089179
364          17.089179
365          17.089179
366          17.089179
367          17.089179
368          17.089179
369          17.089179
370          17.089179
371          17.089179
372          17.089179
373          17.089179
374          19.341350
375          18.164283
376          17.256104
377          17.562815
378          17.562815
379          17.717395
380          16.969341
381          18.312906
382          18.312906
383          18.312906
384          21.535150
385          21.535150
386          23.445404
387          21.897079
388          21.897079
389          20.454256
390          20.454256
391          20.454256
392          20.454256
393          20.454256
394          20.454256
395          20.454256
396          20.454256
397          20.454256
398          21.858624
399          21.858624
400          17.411482
401          17.411482
402          17.411482
403          17.411482
404          17.411482
405          20.535268
406          14.427189
407          15.372445
408          16.297346
409          15.680065
410          19.916058
411          18.570304
412          18.562651
413          22.209059
414          17.983579
415          18.736493
416          17.222956
417          19.348183
418          19.348183
419          19.717872
420          22.913727
421          21.743288
422          21.743288
423          22.913727
424          22.913727
425          30.537281
426          30.537281
427          21.813852
428          21.813852
429          17.801528
430          27.348777
431          27.348777
432          40.512566
433          40.512566
434          40.512566
435          40.512566
436          40.512566
437          40.512566
438          40.512566
439          40.512566
440          40.512566
441          40.512566
442          40.512566
443          40.512566
444          40.512566
445          40.512566
446          40.512566
447          40.512566
448          40.512566
449          40.512566
450          40.512566
451          40.512566
452          40.512566
453          40.512566
454          40.512566
455          40.512566
456          40.512566
457          40.512566
458          40.512566
459          40.512566
460          40.512566
461          40.512566
462          40.512566
463          40.512566
464          40.512566
465          40.512566
466          40.512566
467          40.512566
468          40.512566
469          40.512566
470          40.512566
471          40.512566
472          40.512566
473          40.512566
474          40.512566
475          40.512566
476          40.512566
477          40.512566
478          40.512566
479          40.512566
480          40.512566
481          40.512566
482          40.512566
483          40.512566
484          40.512566
485          40.512566
486          40.512566
487          40.512566
488          40.512566
489          40.512566
490          40.512566
491          40.512566
492          40.512566
493          40.512566
494          20.746522
495          16.288450
496          20.663572
497          20.663572
498          20.663572
499          20.663572
500          20.663572
501          20.663572
502          20.663572
503          20.663572
504          20.663572
505          20.663572
506          20.663572
507          20.663572
508          20.663572
509          20.663572
510          20.663572
511          20.663572
512          20.663572
513          20.663572
514          20.663572
515          20.663572
516          20.663572
517          20.663572
518          20.663572
519          20.663572
520          20.663572
521          20.663572
522          20.663572
523          20.663572
524          20.663572
525          20.663572
526          20.663572
527          20.663572
528          20.663572
529          20.663572
530          20.663572
531          20.663572
532          20.663572
533          19.979493
534          19.166006
535          45.182032
536          20.539460
537          17.631342
538          17.631342
539          47.882099
540          47.882099
541          73.589262
542          73.589262
543          73.589262
544          73.589262
545          73.589262
546          73.589262
547          11.840302
548          12.647299
549          12.857020
550          11.647801
551          73.589262
552          10.433070
553          10.433070
554          10.433070
555          11.706439
556          11.768395
557          13.812919
558          13.181219
559          13.181219
560          10.302434
561          10.541067
562          12.768895
563          11.578704
564          12.523621
565          10.051243
566          10.421326
567          17.153326
568          17.153326
569          15.672649
570           9.266242
571          16.646831
572          16.646831
573          16.646831
574          16.646831
575          12.905475
576          12.407543
577          16.960057
578          15.980155
579          17.780769
580          15.583865
581          10.574860
582           7.205204
583           7.148086
584           8.253521
585           9.488469
586           7.547831
587          10.133495
588           9.373530
589           8.371785
590           9.089013
591          12.221164
592           7.702673
593          12.359854
594          12.605871
595           9.151539
596          10.162648
597          10.048579
598           6.830052
599          10.529484
600          13.249704
601          16.461415
602          16.461415
603          16.461415
604          16.461415
605          16.475316
606          14.186260
607          15.810996
608          15.233389
609          16.513781
610          16.038754
611          15.406476
612          16.039744
613          14.295303
614          14.195876
615          15.080113
616          14.757082
617          14.543009
618          11.283165
619          12.876242
620          13.012346
621          15.197850
622          13.364612
623          13.240551
624          13.779121
625          13.181617
626          25.248076
627          12.198407
628          10.112837
629          20.990581
630          20.990581
631          19.565107
632          19.565107
633          19.565107
634          19.565107
635          12.963362
636          10.614761
637          11.102207
638          10.274178
639          14.738221
640          13.982151
641          16.012976
642          13.991561
643          14.256662
644          33.902975
645          33.902975
646          26.667013
647          29.907590
648          30.293411
649          29.515364
650          31.935244
651          33.816946
652          22.792076
653          29.907590
654          24.850357
655          28.317023
656          29.778303
657          28.772787
658          23.612481
659          26.813957
660          28.542999
661          24.654967
662          24.943135
663          14.356966
664          14.356966
665          14.252179
666          10.961119
667          16.305189
668          15.025566
669          11.938093
670          17.232891
671          12.909037
672          11.427570
673          14.348731
674          12.388851
675          12.388851
676          12.388851
677          12.388851
678          12.388851
679          12.388851
680          12.388851
681          12.388851
682          12.388851
683          12.388851
684          11.972094
685          33.254181
686          33.254181
687          33.254181
688          33.254181
689          33.936259
690          33.492496
691          34.626248
692          29.369544
693          26.187613
694          32.468107
695          34.315131
696          32.106740
697          25.710142
698          13.913381
699          21.368421
700          15.193740
701          14.551417
702          29.867610
703          30.707543
704          29.293516
705          20.300954
706          15.193740
707          12.047799
708          13.548983
709          13.548983
710          25.744869
711          25.568477
712          30.504252
713          13.099185
714          21.090701
715          21.660643
716          12.099604
717          12.099604
718          12.278068
719          18.819584
720          25.166309
721          19.495035
722          19.939265
723          19.777921
724           9.506197
725          14.234373
726          11.298847
727          20.089796
728          11.190665
729          18.165002
730          17.678434
731          17.678434
732          17.388715
733          15.187220
734          15.050692
735          15.183775
736          14.281377
737           9.193501
738          13.284951
739          13.284951
740          23.810077
741          14.283475
742          14.283475
743          14.283475
744          14.283475
745          14.283475
746          14.283475
747          14.283475
748          14.283475
749          14.283475
750          14.283475
751          13.645311
752          35.620812
753          34.177627
754          20.633944
755          20.633944
756          20.633944
757          20.633944
758          36.284382
759          11.822191
760          11.822191
761          11.822191
762          11.822191
763          11.822191
764          11.822191
765          11.822191
766          11.822191
767          11.822191
768          11.822191
769          14.166298
770          37.115128
771          13.457708
772          23.590543
773          34.346713
774          32.700959
775          38.021027
776          38.021027
777          38.021027
778          38.021027
779          38.021027
780          38.021027
781          38.021027
782          38.021027
783          43.992838
784          43.992838
785          18.631403
786          18.430563
787          18.213301
788          18.430563
789          18.430563
790          19.101370
791          18.708063
792          18.631403
793          19.400540
794          19.155325
795          19.011379
796          19.011379
797          26.337160
798          20.580106
799          20.580106
800          20.453848
801          35.110913
802          35.017901
803          34.731418
804          34.731418
805          34.731418
806          34.429487
807          12.305978
808          12.305978
809          12.305978
810          12.305978
811          12.305978
812          12.305978
813          37.965491
814          37.965491
815          37.965491
816          19.970065
817          21.782446
818          17.565326
819          17.759829
820          12.979011
821          13.631427
822          26.179968
823          15.570621
824           9.735484
825          15.765710
826          14.742520
827          16.302873
828          14.325744
829          13.597943
830          15.662600
831          15.662600
832          15.662600
833          14.273696
834          16.615454
835          16.615454
836          18.139917
837          20.007072
838          20.007072
839          19.866546
840          16.338006
841          20.138700
842          23.486569
843          20.238439
844          22.517622
845          22.337084
846          21.817199
847          23.740825
848          21.872604
849          18.735760
850          16.698366
851          16.702122
852          12.071726
853          16.784577
854          14.664233
855          12.866827
856          12.830051
857          19.074949
858          22.005967
859          19.436823
860          19.436823
861          19.588928
862          19.588928
863          17.088811
864          15.163353
865          18.137910
866          21.761596
867          19.450480
868          19.683171
869          18.026928
870          18.283826
871          20.032893
872          17.340288
873          20.124846
874          21.680103
875          26.287342
876          18.997524

Raster aggregation

We can use the same procedure to aggregate a raster dataset to a vector polygon dataset. That’s a very common use case. Let’s reload our Cologne shapefile.

cologne <- 
  sf::st_read(
    "./data/cologne.shp",
    quiet = TRUE
  ) |> 
  sf::st_transform(25832)

plot(cologne["NAME"])

Add the aggregated data

cologne$immigrant_rate <-
  terra::extract(
    immigrant_rate, 
    cologne, 
    fun = mean, 
    na.rm = TRUE, 
    ID = FALSE
  ) |> 
  unlist()

plot(cologne["immigrant_rate"])

5. Conversion & Analysis

Raster to points

raster_now_points <-
  immigrant_rate |> 
  terra::as.points()

plot(raster_now_points)

Points to raster

raster_target_layer <- 
  terra::ext(raster_now_points) |> 
  terra::rast(res = 100)

points_now_raster <- 
  raster_now_points |> 
  terra::rasterize(
    raster_target_layer, 
    field = "immigrants_cologne", 
    fun = "mean",
    background = 0
  )

plot(points_now_raster)

Raster to polygons

polygon_raster <-
  immigrant_rate |>  
  terra::as.polygons() |> 
  sf::st_as_sf()

plot(polygon_raster)

Analysis application: creating a quick ‘heatmap’

OpenStreetMap points data are nice for analyzing urban infrastructure. Let’s draw a quick ‘heatmap’ using obersvation densities.

echarging_sf_cologne_raster <-
  terra::rast(
    echarging_sf_cologne, 
    res = 1000
  )

echarging_sf_cologne_densities <- 
  echarging_sf_cologne |> 
  terra::rasterize(
    echarging_sf_cologne_raster, 
    fun = length, 
    background = 0
  )

plot(echarging_sf_cologne_densities)

Backup / For Home Use

Focal statistics / spatial filter

Focal statistics are another method of including information near a point in space. However, it’s applied to the whole dataset and is independent of arbitrary points we project onto a map.

  • relates focal cells to surrounding cells
  • vastly used in image processing
  • but also applicable in social science research, as we will see

Analysis: edge detection

Source: https://en.wikipedia.org/wiki/Sobel_operator

Edges of immigrant rates

We can do that as well using a sobel filter

\[r_x = \begin{bmatrix}1 & 0 & -1 \\2 & 0 & -2 \\1 & 0 & -1\end{bmatrix} \times raster\_file \\r_y = \begin{bmatrix}1 & 2 & 1 \\0 & 0 & 0 \\-1 & -2 & -1\end{bmatrix}\times raster\_file \\r_{xy} = \sqrt{r_{x}^2 + r_{y}^2}\]

Implementation in R

From the official documentation:

sobel <- function(r) {
  fy <- matrix(c(1, 0, -1, 2, 0, -2, 1, 0, -1), nrow = 3)
  fx <- matrix(c(-1, -2, -1, 0, 0, 0, 1, 2, 1) , nrow = 3)
  rx <- terra::focal(r, fx)
  ry <- terra::focal(r, fy)
  sqrt(rx^2 + ry^2)
}

Data preparation and application of filter

old_extent <- terra::ext(immigrant_rate) 
new_extent <- old_extent + c(10000, -10000, 10000, -10000)

smaller_immigrant_rate <-
  immigrant_rate |> 
  terra::crop(new_extent)

smaller_immigrant_rate[smaller_immigrant_rate < 10] <- NA

immigrant_edges <- sobel(smaller_immigrant_rate)

Comparison

plot(smaller_immigrant_rate)
plot(immigrant_edges)